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Abstract 


This  study  develops  an  analytic  function  for  the  six  dimensional 
surface  (or  hypersurface)  above  the  planet  Venus  for  a  five  earth  year 
survivability  for  an  artificial  satellite. 

Current  US  policy  concerning  the  exploration  of  other  planets,  via 
artificial  satellites,  requires  the  satellites  be  sterilized  (5:61).  This  is  a 
very  time  intensive  and  costly  practice.  Developing  the  ability  to 
estimate  the  life  time  of  an  artificial  satellite  that  can  no  longer  perform 
its  station  keeping  duties  may  allow  the  sterilization  procedures  before 
launch  to  be  waived.  The  objective  is  to  develop  a  five  year 
survivability  function  (denoted  by  hp)  in  the  orbital  parameter  space 
above  which  a  satellite  has  at  least  five  years  to  survive  before  it 
impacts  the  planet’s  surface. 

Perturbations  effects  which  would  cause  the  satellite’s  orbit  to 
deteriorate  are  modeled  and  include:  a)  the  geopotential  of  the  planet,  b) 
the  effects  of  solar  wind,  and  c)  the  drag  on  the  satellite  due  to  the 
atmosphere  of  the  planet. 

The  model  was  then  interpolated  to  provide  an  analytic  function  for 
five  year  survivability. 
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HYPERSURFACE  INSERTION  WINDOW  FOR  LONG 


TERM  ORBITAL  STABILITY  OF  ARTIFICIAL 
SATELLITES  ABOUT  THE  PLANET  VENUS 

L  Introduction 

Currently,  the  United  States  spends  millions  of  dollars  per  satellite 
for  sterilization  prior  to  launch.  There  is  a  concern  that,  should  an 
earth  borne  satellite  impact  another  planet,  the  satellite  would 
contaminate  the  indigenous  environment.  A  likely  scenario  is  a  satellite 
in  planetary  orbit  with  its  station  keeping  fuel  depleted.  Predicting  the 
remaining  life  time  of  the  satellite  (i.e.,  the  amount  of  time  prior  to  the 
satellite  entering  the  planet’s  atmosphere)  could  eliminate  the  need  for 
sterilization  and  replace  this  practice  with  a  rescue  mission. 

This  thesis  predicts  the  life  time  of  a  nonmaneuverable  satellite  by 
generating  a  five  year  hypersurface.  This  is  a  six  dimensional  surface, 
in  the  orbital  parameters  space,  which  marks  the  point  of  five  (earth) 
years  to  impact.  A  satellite  above  this  surface  has  five  years  to  survive. 
Once  a  satellite  contacts  this  surface  there  will  be  a  specific  amount  of 
time  to  initiate  a  rescue  (or  destroy)  mission. 

Specifically,  a  perturbation  model  for  a  satellite  orbiting  Venus  will 
be  calculated  using  a  4x4  geopotential  model  for  Venus,  a  solar  wind 
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model,  and  a  model  for  atmospheric  drag.  In  order  to  decrease  the  size 
of  the  dimensional  space  the  initial  mean  anomaly  M  and  the  initial 
argument  of  periapsis  w  will  be  considered  constant.  This  will  reduce 
the  space  to  four  dimensions. 
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IL  Model 

In  generating  the  model  for  this  thesis,  there  were  three 
perturbations  considered:  the  geopotential  of  Venus,  atmospheric  drag, 
and  solar  wind.  The  perturbations  are  given  in  the  following  sections 
along  with  Lagrange’s  Planetary  Equations  (LPE)  and  a  model  of  a 
typical  satellite. 

Lagrange’s  equations  and  the  perturbations  are  given  in  terms  of  a 
hybrid  set  of  parameters.  For  low  eccentricities  (on  the  order  of  0.02  to 
0.10)  the  eccentricity  and  the  argument  of  perigee  are  replaced  by  two 
equinoctial  elements,  h  and  k  (from  a  coordinate  system  with 
singularities  at  i  -  n  and  for  rectilinear  orbits).  Also,  for  long  term 
effects  of  the  orbit  the  short  period  mean  anomaly  M  is  replaced  with 
the  stroboscopic  mean  node  (A*).  This  stroboscopic  mean  node  acts  as  a 
discrete  measurement  of  the  variation  of  the  mean  anomaly  at  times  df. 
The  inclination,  the  longitude  of  the  ascending  node,  and  the  semimajor 
axis  are  given  typical  values  for  an  experimental  satellite.  This  set  of 
orbital  parameters,  along  with  their  relationship  to  the  classical  set,  will 
be  explained  in  more  detail  in  the  next  section. 

Lagrange’s  Orbital  Equations 

Lagrange’s  Planetary  Equations  (LPE),  written  in  terms  of  the 
orbital  elements  (a,e,i,ft/o,M),  are  (13:29): 
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da  2  dR' 
dt  nadhd 

do  1  -  o2dR‘  Jl-o2dR‘ 
dt  na2od A/  na2o  daj 


dev  cost  dR‘  V  1  -  o 2 d R ' 

- m  —  - -■  — - - ♦ - 

dt  na2J\-  e2sin i  ^  na2e  do 

di  cost  dR’  1  dR’ 

dt  na2J  1  -  e:sin i  na2J  l  -  o2sini  ^ 


d*  na2J l  -  c2sin  i 


dW  l-e2dK*  2  a** 

— —  -  /i - - — - - — 

dt  na2 a  3o  na  da 


(2.1.1a) 


(2.1.16) 


(2.1.1c) 


(2. 1 .Id) 


(2.1.1c) 


(2.1.1/) 


Where  /?*is  the  disturbing  potential  and  n  is  the  mean  motion. 


The  LPE  contain  terms  with  the  eccentricity  and  the  inclination  in 
the  denominator.  For  the  small  values  of  these  elements  considered  here 
(0.02  <«  SO.  10  and  0.5  S<S  15.5)  a  more  well  behaved  set  of  variables  will 
be  used:  two  of  the  equinoctial  elements,  h  and  k.  These  elements 
eliminate  the  singularity  due  to  zero  eccentricity,  and  have  their  own 
singularity  away  from  the  range  of  interest  at  i  -  n  (22:23).  h  and  k  are 
given  by 


h  -  esin  o> 


(2.1 .2a) 
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k  =  B  COS  60 


(2.1  .26) 


Since  only  long  term  orbital  behavior  is  of  interest,  the  fast 
variable  M  (the  Mean  Anomaly)  will  be  replaced  with  the  stroboscopic 
mean  node,  A»(7:167-189).  Resonance  occurs  when  integer  multiples  of  the 
mean  secular  rate  of  the  various  orbital  parameters  are  matched  with 
the  coefficients  of  the  geopotential.  These  resonance  effects  are  due  to 
longitude  dependent  tesserals  in  the  central  body  gravity  field.  The 
stroboscopic  mean  node  reveals  the  most  pronounced  effects  of 
resonance.  This  allows  us  to  retain  any  resonance  effects  due  to  tesseral 
harmonics  (13,1966:49-56).  The  stroboscopic  mean  node  is  given  by 

(2.1.3) 

where  50  is  the  ratio  of  two  relatively  prime  integers,  P/Q, 
approximating  the  number  of  nodal  crossings  per  planet  revolution.  9,  is 
the  prime  meridian  angle.  Gedeon  (7:171)  defines  A*  by  considering  a 
"mean  satellite".  Now,  input  at  fc=0.  If  Venus  is  illuminated  with  a 

strobe  light  the  mean  satellite  will  be  seen  above  the  equator  at 
A „  =  n-8,  longitude.  If  Venus  is  flashed  again  after  Q  days,  then  (if  A*,*0 
)  the  mean  satellite  will  be  at  the  same  longitude.  However,  if  A„*0,  the 
satellite  will  be  a  distance  /A„d(  away.  This  gives  a  measure  of  the 
change  of  the  mean  node  from  period  to  period. 

The  LPE  must  be  written  in  terms  of  the  new  orbital  parameter  set 
(a,h,i,k ).  First  transform  the  disturbing  potential  /?( a,e,i,A,o>,M)  into 
a  function  of  the  new  elements  A(a,h,i,kAA«). 
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dR‘  LdR  dR  1  dR 

- •  k - h  —  + - 

daj  dh  dk  S0dkN 

dR'  _  hdjl  kdR. 

do  o  dh  o  dk 


dR *  dR  dR 
dn  "  dCl*  dk~N 

dR'  1  dR 
dM  ™  Sodkf, 


(2.1  .4a) 


(2.1  .46) 


(2.1.4c) 


(2.1  .4d) 


The  orbital  equations  for  h,  k,  and  are 


dh  hdo  dco 
dt  o  dt  dt 


(2.1  .5a) 


dk^kdo  d  go 
dt  odt  dt 


(2.1  .56) 


dk,, _  l  f dM  +d(v\dfl  dOt 
~dfm  Jl\~dT*  ~dT  )  ~dt~  ~dt 


(2.1 .5c) 


Now  rewrite  the  LPE  into  the  new  set  of  orbital  parameters  by  way 
of  the  above  equations(9:2)  (15:1). 


da  2  1  dR 

—  —  a - ...  — . 

dt  naSoJA# 


(2.1.6a) 


dh_yl\-o*dR  it  cot  t  dR  hJ\-o2  dR 
dt  na 2  dk  na*J  \  -  di  na2S0 


(2.1  .66) 


6 


where  o 


1 

1  ♦  V 1  ~  e2 


Potential  MsdsI 


In  order  to  get  an  expression  for  the  gravitational  potential 
(geopotential),  first  look  at  the  Laplacian  of  the  potential.  Assume  the 
geopotential  to  be  conservative,  then  the  Laplacian  will  be  zero.  In 
tensor  notation 
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where  gkl  -  contravariant  metric  tensor 
/g-  [det  (g„)],/a 

Using  spherical  coordinates  (radius (r),  longitude(t),  and  latitude(ft)), 
define  the  differential  length  as 

(ds)2  •  (dr)2  *  (rcosBdj)2*  ( rdO )2 
This  will  make  the  covariant  metric  tensor 

n  o  o\ 

[g(/]-  0  r3cos20  0  I 

VO  0  r2J 

and  the  contravariant  metric  tensor 


Now,  to  evaluate  the  Laplacian,  note 

/g  -  r2cos0 

and 

i#  j  +  gl‘m  0 
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For  i=7=l:  x1  =  r  and 

1  a 


a  Jv 
r  cos0  — 
br 


r2 cos03r 

For  i=j=2:  x*  =  4  and 


1.1±( 

r2arl 

>  dr) 

I  a  f  ,  „  1  bVl  1  a2i 

— - r  COS  0  — r - — -  I  ■  — - 

r  cosea^L  r2cosfla0j  r  cosfli^ 


For  i'=j'=3:  x3  =  0  and 

1  a 


r3cos0b0 


2  a  1  *V 
r  COS  0  r  — — 
r2b0 


1  a 

r3cos8b9 


(cos9j?) 


Substituting  eqn.s(2.2.2)  into  the  Laplacian  of  the 
(eqn(2.2.1))  and  assuming  a  conservative  field  yields 


1  b3V  1  a  r  av\ 

r2cos00^2  r3cos$b9\  b8 ) 


7V-0 

Assuming  the  potential  is  linear  and  separable,  then 
K(r,0,*)-*(r)*(9)*(*) 


and  substituting  eqn(2.2.4)  into  eqn(2.2.3)  yields 


(2.2.2a) 

(2.2.2b) 

(2.2.2c) 

geopotential 

(2.2.3) 

(2.2.4) 
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♦  R* 


1  d  (  ndB\ 
- cos  0 - 

cos0d0V  d6  J 


*  RB 


( 


1  d24> 

cos  20dp2 


) 


0  (2.2.5) 


Multiply  the  above  equation  by  cos2  9/RBt  to  get  Laplace’s  equation 
in  the  form 


cos28  d  f  2dR  \  ^ 
R  dr\  dr  j 


cos  0  d  ( 

~lTde{ 


COS0 


dB\ 

dej 


+ 


1  d2<t> 
♦  di>2 


(2.2.6) 


Since  *  is  alone  in  Laplace's  equation,  it  can  be  separated  from  the 
rest  of  the  equation  and  both  sides  can  be  set  equal  to  the  separation 
constant  -m2,  so 


d2» 

d*‘ 


■+  m  0 


(2.2.7) 


This  is  a  Sturm-Liouville  equation  with  the  solution 


$($)■  Cmcosm4>+  Smsinm4>  where  m-  1,2,3,. 


(2.2.8) 


Rewriting  Laplace's  equation 


cos20 


f  d  (dR\~\  cos 0 {"  d  ( 

LdrldrJj+  B  L  0  V  C°S  ^  d  0  j  J  " 


m 


and  dividing  by  cos2  8  and  separating,  produces 


Rdr\  dr) 


m 


1  d 


cos20  BcoeBdO 


(2.2.9) 
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To  evaluate  Bfi),  using  eqn(2.2.9),  first  rewrite  the  right  hand  side 


1  d 
cos  Odd 


„dB\  (  m 2  \ 

id—  I-  I  - — -  /)  fi  -  0 

dd  J  \cos30  ) 


Next  transform  the  above  equation  to  the  form  of  a  Legendre's 
Associated  ordinary  differential  equation.  Let  sin®  =  x  and  1  *  n(n+l)  so 


cosOdd 


^(cos95fHn(n'nTS]fl‘° 


(2.2.10) 


where  the  identity  cos2®  =  1  -  sin2®  =  1  -  x2  was  used. 


Evaluating  the  derivatives  yields 


1  d  f  dB 

cosedel  de 


)  -  dB  ,,  sd!* 


Substituting  into  the  Associated  Legendre's  equation,  eqn(2.2.10): 


.  2\d2 B  „  dB  f  ,  lx  m2  1  „ 

1  -  x2 — r-  2x  —  +  n(n  +  1) - -  \B  ■  0 

'  1 dx 2  dx  L  1  - x 2 J 


(2.2.11) 


The  solution  to  the  Associated  Legendre's  equation  is 


B(0)m  /'"(sin 0) 


(2.2.12a) 


where  l  -  n(n+  1 ) 


P?(  8in0)-cos"0  £  Tlmt8in‘"n'2t0 

t^o 


(2.2.12b) 
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T  (-1),(2/-2Q! 

2'/ !(/-/)!  (/-m- 20* 


(2.2.12c) 


All  that  is  left  is  to  solve  for  R(r).  Looking  at  the  potential 
expression  where  B(0)  was  separated  (eqn(2.2.9))  yields 

r2^r)-n(n*  \)R-0  (2.2.13) 

dr  v.  dr  J 


There  are  two  possible  solutions  to  the  above  equation.  The  first 
solution,  R(r)  =  rn,  blows  up  as  r->  ®.  Therefore,  only  keep  the  solution 

R(r)-rl'“,)  (2.2.14) 

Combining  the  solutions  for  R,  B,  and  ♦  generates  the  following 
relationship  for  the  geopotential  (in  terms  of  radius,  longitude,  and 
latitude).  Note  that  the  terms  m  and  R»  have  been  added.  This  is  to  make 
the  Cim  and  Sjm  terms  dimensionless. 

V(r  .0.*)- ?r(sin0)(C  ,mcosm^«-  S  ,„sinm^)  (2.2.15) 

r  >-0m?0\  "tj 

Now  translate  the  geopotential  from  its  spherical  harmonic 
representation  (eqn(2.2.15))  to  a  Keplerian  representation  (see  Appendix 
A  for  a  full  derivation).  In  order  to  do  this,  two  new  functions  must  be 
introduced:  the  inclination  Function  Fimp(i )  and  the  Eccentricity  Function 
G^q(e). 

The  Inclination  Function  is  defined  as 
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<•>... 


_  ...  c-  zi 

*  inf  (0-  )  1 1  / 1  77771 


(  2  /  —  2  f )  •  t>< 


2t)!  221'2' 


sin11  "  "’(i) 


The  Eccentricity  Function  may  be  defined  as 


£  tf  (ip-n  " 


a  '*' ( 1  -  e*) 


*-2p']l  d  Jl  2  J  (2>2'  ) 


where 


' p  for  pi  1/2 
J-  p  for  pZ  1/2 


Hence,  the  disturbing  function  for  a  nonspherical  planet  is  given 


by 


•  l 


1-2  m-0 


(2.2. 18a) 


where  V ,m  -  ^  £  Flm,(i)  £  G,Pf(o)5 ,npt((o  ,hi  ,  fl,0) 

a  p-o 


**-« 


(2.2.18b) 


I C  lnco8j+  S  jmsin^  ,  f-n 

■■Id 


-  m  even 
m  odd 


(2.2.18c) 


where  ^  =  (1  -  2p)ou  +  (/ - 2p  ♦  q)M  *  m(fi  -  9)  (2.2. 18d) 


Or,  using  the  conversions  of  Appendix  B,  write  V*  as 
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( 5.6a ) 


tm  1 


i  - 


!•  I  J  ,m 


X  5  <mp(0  X  ^  ’*1 

P  *0  *  —  - 


») 


where 


<5  (mp  8  °  Imf 

.  1  (cosf  eU!cosqco  +  sin f  o'*'sin geo  ,  l- m  even 

ewl \sin  £  eu,cosga>  -  cosf  e'^sin  qoo  ,  l  -  m  odd 


4>‘  -  (l  -  2  p)u)  +  (l  -  2 p  +  q)M  +  m((l-  Q)~  mk 
Calculating  the  derivatives  of  the  geopotential  yields 


&V  im  (/+l)/i5ir  V  r  r*  c 

aa  a'*2  J,mh'  ,nf 


Ipt 


{C- 

P-0 


dh. 

a'*1 

'N 

3 

1 

di 

a 

*V,mm 

a*5: 

dk 

a'*1 

i 

6 

pR\ 

an 

a'*1 

(PP  a/l 


hdC') 
e  de  f 


l  t 


dF 


Imp 


£  ipp^ipf 


dS 


IWif 


_ a  C 

(p,  dk 


kdG'\ 
o  do  j 


PR  \  ,  X~  Imp  n  .  1  25  Imi  _  ^  A 


(5.66) 

(5.3) 

(5.4) 
(5.1c) 

(2.2. 19a) 

(2.2.196) 

(2.2.19c) 

(2.2. 19d) 

(2.2.19c) 
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bV 


Im 


bn 


0 ,  m 


0 


bVlin 


-  0 


These  results  (eqn.s(2.2.19))  may  be  inserted  into 
Planetary  Equations  (eqn.s(2.1.6)l  resulting  in 


(2.2.19 /) 

(2.2. 1 9  gr ) 

Lagrange’s 

(2.2.20a) 

(2.2.20b) 

(2.2.20c) 


15 


hcoti  pR1,  ,  o 

*  2  f-  j  I-. L  1,  ^ 

na  v  1  ■  «2q  *-o 


(2 . 2 . 20d) 


dt  m  2  r. - i  C  (2.2.20e) 

na  vl“®  sinia  *-o  £,i 


dA.*,  n  d0,  1 

•  + 


df  5 o  d/  na2S 


-[Vi-«2/r 


■  (***•  i  V  / 1,  r  f,~*  *Slinl.e>  hdC'} 

U'*'  'mp-c\ ""'I  *M  ^  ‘"'a  J 


+  *f 


Imp 


(c-„. 


a*  3"n*cdo,/) 


a  /?  o  ~  cos  I  fi  r.  ,  _  r-  a  f  i 


'2q^‘  77'  ■■■;-  -^»Eirc'Wf'-  I  (2.2.20f) 

<>a  ^  1  -  a2sin  i°  »-o 


] 


Solar  Wind  Model 

The  acceleration  of  a  satellite  due  to  direct  solar  radiation  pressure 
is  given  by(8:12) 

(2-31> 

where  r-  position  vector  of  satellite  relative  to  planet  center 
Ct««*  distance  of  sun  at  1  A.U. 
rs  -  position  of  sun  relative  to  planet  center 


16 


r„«  r- r5 

pt-  radiation  pressure  on  a  perfectly  absorbing  surface 
y-  black  body  radiation  constant  (-  1) 

m  =  mass  of  the  satellite 

A'  =  satellite  reference  area 
r  =  unit  vector  designation 

At  times  Venus  will  be  between  the  sun  and  the  satellite  with  Venus 
blocking  the  solar  radiation  and,  hence,  the  need  to  compute  the  region 
for  this  solar  occultation  g(u)  is  necessary.  To  do  this,  a  right  handed 
coordinate  system  (Si,S2,S3)  will  be  employed  whose  origin  is  at  the 
center  of  Venus.  Si  will  point  at  the  sun  and  S3  will  point  "up"  out  of 
the  plane  of  the  ecliptic.  With  this  frame  of  reference,  the  conditions  for 
occultation  are 


5,<0  (2.3.2a) 

and 

S22+S;<R2f  (2.3.2b) 

where  Rp  =  planet  radius  plus  atmospheric  altitude  blockage.  If 

r2  ■  +•  5“, +  5“, 

then  eqn(2.3.2b)  may  be  rewritten  as 
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r2-S*<R2t 


or 


Syl-JP-Rl  (2.3.3) 

(i.e.,  Si  must  be  more  negative).  Now  eqn(2.3.3)  may  replace  eqn.s(2.3.2) 
as  the  occultation  condition. 

Next,  assume  the  right  ascension  («)  and  declination  (6)  are  known 
in  some  reference  frame.  Then,  the  cartesian  coordinates  of  the  satellite 
are  given  by 


x  i  *  r(cosjQcosu  -  cos i sin  (1  sin u) 


x  2  *  r(sinf2cosu  ♦  cosicosHsin  u) 


xz  •  r  sin  isin  u 


(2.3.4a) 

(2.3.4b) 

(2.3.4c) 


where  In  order  to  relate  the  Si,Si,S3-coordinate  system  to  the 

cartesian  coordinates  of  the  satellite,  two  rotations  will  be  performed: 


cos(-6)  0  -sin(-6)' 

0  1  0 
sin(-6)  0  cos(-6) 


cosa  sina  ON  f  x , 

x|  -  sin  a  cosa  0  x2 

0  0  1  J\x3j 


or 
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cosficosa 

cosdsin  a 

sin 

r 

-  sin  a 

cosa 

° 

* 

-  sin  6cosa 

-  sin 6sin a 

cos6/ 

i* 

Using  eqn(2.3.5),  eqn(2.3.3)  can  be  written  as 

Si  »  (cosbcosa)x  |  +  (cos6sina)x2  ♦  (sin  6)x, 

For  ease  of  manipulation,  define 

A  •  cosfisina 

B  •  sin  6 

C  •  cosficosa 

Using  eqn.s(2.3.4)  and  eqn.s(2.3.7),  eqn(2.3.6)  becomes 
S ,  -  rC( cos fl cos u  -  cos  i  sin  (1  sin  u) 

♦  r^(sin  Clcosu  *  cosicos/2sin  u) 

♦  r£(  sin  isin  u) 

Again  for  simplicity,  define 
F  •  C  cos  Cl*  A  sin  Cl 
C  •  cosiMcosfl-  Csinfl)*  Bsini 


(2.3.5) 


Rl  (2.3.6) 

(2.3.7a) 

(2.3.7b) 

(2.3.7c) 


(2.3.8) 

(2.3.9a) 

(2.3.9b) 
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H  - 


(2.3.9c) 


Using  eqn.s(2.3.9),  the  occultation  condition  is 

g(  u )  =  F  cosu  *  G  sin  u  *  //(a)  <  0  (2.3.10) 

There  are  three  restrictions  to  the  above  occultation  condition, 
gfu).  First,  g(u)  depends  on  cosu  and  sinu  causing  it  to  be  periodic 
with  period  2n.  Second,  g(u)  also  depends  on  H(u)  implying  the  minimum 
of  g(u)  may  not  be  zero.  Finally,  when  gfu)  is  at  a  minimum  and  u=umin 
then  the  satellite’s  orbit  is  occulted  at  least  part  of  the  time  from  the 
sun. 


The  following  is  necessary  in  order  to  calculate  the  orbital 
equations  due  to  6olar  radiation.  Since  the  orbital  eccentricity  of  Venus 
is  on  the  order  of  0.0068,  it  can  be  assumed  i*  =  ap  (where 
ap  =  semi-major  axis  of  the  Venus  orbit).  Also  assume  ivun  =  at 
(a>  =  semi-major  axis  of  the  earth  orbit,  e*  =  0.0167).  Since  the  satellite 
is  in  orbit,  assume  r,"-r,  and  the  radiation  pressure  equation  becomes 


•  C'A' 

c - 

m 


(2.3.11) 


where  C, 


yp 


2 


Note  that  the  above  equation  is  of  the  form  of  a  gradient  of  a 
potential  (R): 
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r-VR' 


where  /?’ 


CtA' 


m 


C 


or  R 


Ct  A ' 


m 


rcosS 


(2.3.12) 


Where  5  =  planet  centered  angle  between  the  sun  and  the  satellite.  The 
cos S  may  be  determined  (from  eqn(2.3.8})  to  be 


cos  S  -  C(cos/2  ♦  cosu  -  cos  isin  /Isin  u) 

+  j4(sin /lcosu  +  cost  cos /2  sin  u) 

♦B(sintsinu) 

The  potential,  R  *,  is  in  the  classical  orbital  set  and  has  the 
following  derivatives 


dR’  C'A'r 

— — - - -  —  cos  S 

da  m  a 

dR ’  CtA'  dr  C,A'  d(cosS)3u 

- - cos  S- - r — - - — 

da  m  do  m  du  do 

dR *  CtA  T  dr  d(cosS)du~ 

— - cos  5  —  +  r — - — 

do  m  [  do  du  do  _ 

dR *  CtA'  d(cosS) 

- - - r - 

di  m  di 


(2.3.13a) 


(2.3. 13b) 


(2.3.13c) 
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bR' 

C(A'  b(cosS) 

bCl  " 

m  '  bCl 

bR'  m 

C,A'  b(cosS) 

b(JD 

m  '  bu 

bR‘  mdR bf 
bM  ~  bf  bU 


bR' 

bM 


C  tA‘ 
m 


. br 

cos  S  —  ♦  r 
bu 


b(cosS) 

bu 


In  order  to  insert  the  above  relationships  into 
equations,  the  following  identities  are  needed: 

br_  r  2csin/ 
bu  a(  1  -  e2) 

br  , 

—  ■  -  a  cos  / 
bo 


it 

bo 


sin 


/ 


3(cos5) 

__ 


C(sin  t  sin  Cl  sin  u) 


♦  ^(  -  sin < cos  A sinu)  +  B (cost sin  u) 


£(cos5) 

bCl 


C(-  cos cos u-  cos  j  cos  Cl  sin  u) 


(2.3. 13d) 

(2.3.13*) 

(2.3. 13f ) 

Lagrange’s  orbital 

(2.3.14a) 

(2.3. I4b) 

(2.3.14c) 

(2.3. 14d) 
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A{-  cos  flu  -  cosj  sin  Cl  sin  u) 


(2.3. 14e) 


d(cosS) 

du 


C(-cos/2sinu-  cos i sin  (1  cos u) 


♦  ^(-sin/2sinu  +  cos  i  cos  12  cos  u)  ♦  B  (sin  icosu)  (2.3. 1 4f ) 

In  the  above  equations,  the  following  partial  derivatives  have  been 
used  interchangeably  (see  eqn(2.4.19)): 

—  m  i  -  6 

df  dcu  du 


Lagrange’s  Planetary  Equations,  in  the  new  variables,  become 


da 

~dt 


—  *dr  .  _^(cos5)' 

cos  o  ~ —  ♦  r  — - - 

du  du 


(2.3. 15a) 


dh 

dt 


na 


r--.-(sinu(2  +  kcosu  *  hsin  u)  +  h)r  *--°-S 
Wl-83 


pj\~02 

- cos  u  cos  S  — 

na 


pkcoti  d(cosS) 
na2J  1  -  e 2 


(2.3. 15b) 


d«  p  cot  i  d(cosS)  p  3(cos5“) 

dt  na2yj  i  -  e2  na2J\~-o2  sini 


dk 

dt 


(cosu(2  +  *  cos u  +  h  sin  u)*  k)r 


d(cosS) 

du 
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p V  I  -  o 2  „  phcoti  £(cos S) 

+ - cosucosS- - ,  - r - — -  (2.3.  ISd) 


na 


na2'j  \  -  q2 


dl 


dd  p  .  .  3(cos5) 

=  sinir 


dt  na 2V 1  ~  e2 


di 


(2.3. 15e) 


dkN  pfi'  d(cosS) 

:(2  +  Jfccosu+/i  sin  u)(k  sinu  -  h.  cos  it  )r 


d*  S0na2yj  1  -  e 2 


(J  u 


S0na 


— — ^(a/?'/T~-T2(*cosu+  hsinu)  +  2r)cos5 

r>  n  -  v  / 


50na 


p  f  S q~  cosi\ 

2V 1  -«H  sin*  J 


where 


C5^' 


m 


Atmospheric  Drag  Model 

The  force  (D)  per  unit  mass  due  to  atmospheric  drag  is  given  by 
(21:81) 

(2.4.1) 
2  nr 

where  0=  drag 

m  -  mass 

C,-  drag  coefficient 
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p«  atmospheric  density  =  p(r) 

A  •  cross  sectional  area  of  the  satellite 

v-  satellite  velocity 

Using  McCuskey’s  notation  (19:81),  the  drag  acts  opposite  to  the 
velocity  of  the  satellite  and  may,  therefore,  be  written  in  vector  notation 
as 


D  - -Dur  (2.4.2) 

where  dr  is  the  unit  vector  along  the  tangent  in  the  direction  of  motion. 

For  ease  of  calculation  it  is  necessary  to  transform  the  drag  in 
terms  of  the  radial  and  transverse  components  of  the  satellite's  flight 
path.  The  radial  unit  vector  (<2,)  is  along  the  position  vector  (r)  of  the 
satellite  orbit.  The  transverse  unit  vector  (fl,)  is  perpendicular  to  c  in 
the  orbital  plane  and  forms  an  acute  angle  with  the  velocity  vector  (t). 
In  this  coordinate  system  the  drag  can  be  written  as 

E-R'ur*S'ut  (2.4.3) 

In  order  to  determine  If  and  S’,  the  following  relationships  are 
employed: 


Dar--D[aT-a,) 

(2.4.4a) 

Dat~-D{aT-at) 

(2.4.4b) 
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To  calculate  the  dot  product  of  the  above  equations  (eqn.s(2.5.4)) 


the  velocity  (v)  may  be  expressed  as 
vm  fur  +  r } u.t  -  vu.T 


so 


r  -  D{a,-G.r)m  -  f-o 


(2.1.5a) 


(2.4.5b) 


S' 


D[&t’  &a)m  ~~~D 


(2. 4. Sc) 


The  following  two  identities(19:148)  will  be  needed: 

1  +  ocos / 

r2f»h-yjk2Ma{\  -  o2) 

Take  the  time  derivative  of  eqn(2.4.6a) 

.  fc/A/  osin  / 

Va(  1  -«J) 

Also,  note(19:148) 

.  ^  A/a7(  1  +  ocos /) 

Va(  1  -«J) 

k<fM(  1  +  ©J*  2 ocos /)l/2 
ly  ■ - ====== - 

V°(i  -•*) 


(2.4.6a) 

(2.4.6b) 


(2.4.7a) 


(2.4.7b) 


(2.4.7c) 
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Divide  eqn(2.4.7a)  by  eqn(2.4.7c)  and  insert  into  eqn(2.4.5b)  to  get 


R' 


-esin  /  D 

( 1  +  a2+  2ecos  f)t/2 


(2.4.8a) 


Divide  eqn(2.4.7b)  by  eqn(2.4.7c)  and  insert  into  eqn(2.4.5c)  to  get 


-  ( 1  +  ocos /)D 
( 1  +  o2  +  2ecos/)l/2 


(2.4.8b) 


Now  place  i?  and  S’  into  the  standard  orbital  reference  frame  with 
a  unit  vector  (P)  along  the  perihelion  and  a  unit  vector  (Q)  at  an  angle 
/  =  90*  to  P.  This  will  make  the  drag 

£-£■/?  'ur  +  S'u, 

or 

F  =  (R'cosf  -  S'sin  f)P+  (R'sin/+  S'cos/)Q  (2.4.9) 

Note  that,  in  this  reference  frame,  the  position  is  given  by 
r  =  rcos  f  P+  r  sin  fQ 


or 


r 


a(l-sJ)cos/-  a(l-oJ)  A 
- p  + - ! - L - Q 

(  1  ♦  ecos /)  ( 1  +  ecos/)sin  / 


(2.4.10) 


Finally,  calculate  the  partial  derivatives  of  the  perturbation  (due  to 
atmospheric  drag)  with  respect  to  the  orbital  elements.  Note 
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bR  _  bR  br, 
be,  bribe, 


or 


bR  r  dr 
be,  ~  be, 


(2.4.11) 


Where  the  cj  are  the  six  classical  orbital  parameters. 


Calculating  the  partial  derivatives  of  the  perturbation: 


~-R'-  (2.4.12a) 

da  a 


b_R 

bo 


-R'a  cos/ *  S  'a  sin/ 


(  1  +  a(  1  -  e2)} 


(2.4. 12b) 


bR  /Tea  sin/  S'a2J  1  -  o2 
bM~  fTTp  +  r 


(2.4.12c) 


bR 

bO 


S'r  cost 


(2.4. 12d) 


bR 

bco 


S'r 


(2.4. 12e) 


bR 

bi 


0 


(2.4. 12f) 
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Insert  eqn.s(2.4.12)  into  Lagrange’s  Planetary  Equations  (eqn(2.1.1))  to 
obtain 


da 

CdA 

dt 

m  P 

do 

CdA 

dt  " 

_  pi 
m 

doo 

- m 

CdA 
- / 

3/2 


1  +  o2  +  2acos  / 

l  J 


.J  l  +  e2+  2ocos/^ 
pna  (o+ cos/ )l - — -  I 


I  /  2 


dt  m 


sin  /  f  1  +  e2  ♦  2ocos/^ 

pnQ—  [ - - J 


I  /2 


(2.4. 13a) 


(2.4.13b) 


(2.4.13c) 


di  . 
—  -  0 
dt 


(2.4. 13d) 


da 

dt 


-0 


(2.4. 13e) 


dM 

dt 


CdA 

m 


pna 


sin  /  ( 1  ♦  o^ocos/) 
o  (1  +  ocos/) 


V  1  +  o22ocos / 


(2.4. 13f) 


To  transform  the  above  equations  into  the  orbital  set  (a,h,i»i<»lV« )  recall 
(from  eqn.s(2.1.5)) 


dh  "hdo  +  ^dou 
dt  o  dt  dt 


(2.4.14a) 


dk  m  kdo  ^  dot) 
dt  q  dt  dt 


(2.4.14b) 
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dkN  _  1  f  dM  +  daj  \  dCl  dOf 
~dt~m  S~0{~dF*  ~dT )  ~dt  ~dt~ 


(2.4.14c) 


Substituting  the  classical  equations  for  atmospheric  drag  (eqn.s(2.4.13)| 
into  the  above  equations  yields 


dh 

dt 


C4A 

m 


pna(o*  cos/) 


^  l  *  o2  +  2ecos/y  2 


+  k 


C4A  sin/ 
■pna - - 


m 


^  1  +  e22ecos  ,J-J 


dh 

dt 


-  4 A  f  1  ♦  e2 *  2ocos/ ^ 

—  - j-p - J 


I  /  2 


(2.4. 15a) 


dk 

dt 


C4A 

m 


pna(o+  cos/) 


C4A 

- pna 

m 


sin  / 

e 


( 

( 


1  »  o2*  2ocos/y/2  j 
1  +  «22ocos/V /2_J 

^  )  J 


dk 

dt 


C4A 

m 


pna(k+ cosu) 


( 


1  ♦  q2+  20C0S/V2 
1-s2  ) 


(2.4.15b) 


d  kN 

hr 


C4A 

- pna 

m 


sin /( 1  +  02  +  geos/) 

0  (  1  ♦  0COS/) 


V  1  +  02cos  / 


C4A 

- pna 

m 


sin  f(  1  +  o2+  2acos/y/2^ 

—[  r^p  )  J 
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dk„ 

dt 


1  CdA 
S0  m 


pnao  sin 


in/  O'---  ■  / 

L  a  J 1  - 


1  +  ea2®cos/  (2.4.15c) 


where  ecos/ *  k  cosu  h  sin  u 
«sin/  =  Asinu-  hcosu 


it  =  to  +  / 


Now,  employ  the  Method  of  Averaging  to  eliminate  the  fast  variables 
and  retain  only  the  averaged  effects  of  atmospheric  drag  over  one 
period.  Introduce  the  averaged  rate,  given  by  (17:11) 


dt  pj  dt 


(2.4.16) 


where  £«a,  h,  k.  or  A,* 

We  can  also  relate  time  to  the  true  anomaly  (/)  by 


dt 


r2df 

na2yj  l  -  q  2 


(2.4.17) 


so  that 

^  n  (2'dK  r2df 

dt  "  2/r  Jo  d<na*/T— ?  (2.4.18) 

The  argument  of  periapsis  (to)  will  have  only  a  small  change  over 
one  period  (see  section  2.5),  so  write 

df^doj  +  du  (2.4.19a) 
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or 


d/«du 


(2.4. 19b) 


Substituting  the  orbital  equations  (eqn.(2.4.13.a)  and  eqn.s(2.4.15)) 
into  the  above  integral  (eqn(2.4.18))  yields 

~m  —  f  pf-1  [  1  +  eJ+  2(fccosu  +  hsinu)V'2du.  (2.4.20a) 
at  i  -  a  Jo  \  a )  i 

dm  r2t  ( rV,. 

7imyh  PUJ  (',*sinu) 


5-0 

dt 


x [  1  +  «2+  2 (* cos u  +  /i sin  u)]l/2du  (2. 4. 20b) 


(2.4.20c) 


s-r-(s)  (*+c»su’ 


x[  1  +  eJ+  2(*cosu  +  hsin  u)]'2du 


dt 


-0 


(2.4.20d) 

(2.4.20e) 


m  y  j  (Jfcsinu-/icosu)[  1  +  a2*  2(jfccosu+  /isinu)] 


1/2 


du 


(2.4.20f ) 
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where  y 


C  t A  na 2  1 
m  1  -  o22n 


p0- reference  density 
h0-  reference  altitude 
h  =  satellite  altitude 

H  =  scale  height 

Resonance  Model 

Whenever  a  parameter  behaves  in  a  sinusoidal  manner  and  is 
effected  by  a  sinusoidal  forcing  function  the  problem  of  resonance  must 
be  addressed.  Resonance  i6  a  phenomenon  whereby  a  displaced 
parameter  approaches  infinity  when  the  frequency  of  the  applied  force 
equals  the  natural  frequency  of  the  parameter  (9:47).  Resonance  effects 
are  potentially  present  in  the  geopotential  model  (13:49)  and  will  be 
described  in  this  section.  This  section  will  also  examine  the  elimination 
of  the  fast  variables  from  the  geopotential. 

The  geopotential  may  be  written  as 

a  p-o  r- 


where 
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,M  .n.o.k 


cost' .  I- rn  even 
sin^’  ,l-m  odd 


(B.lb) 


t’  •  p)oo  *  (/  -  2  p*  q)KI  ♦  m((l-  6)-  mklm  (B  Ac) 

Lagrange's  Planetary  Equations  are  concerned  with  the  time  rate  of 
change  of  the  orbital  elements.  The  geopotential  is  a  forcing  function 
whose  trigonometric  argument  is  given  by  eqn(E.lc),  above. 

Fortunately,  the  rotational  rate  of  Venus  is  retrograde  and  very 
slow  with  a  value  of  -1.71460706e-05*/sec  or  243  days  per  revolution. 
This  means  the  geopotential  will  have  a  frequency  opposite  to  the  orbital 
elements  and  will  be  of  no  concern. 


There  are  three  orbital  elements  in  eqn(B.lc)  which  can  change  with 
time:  M,  u>,  O.  The  mean  anomaly  M  varies  with  the  speed  of  the  orbital 
period  through  the  mean  motion.  A  typical  period  for  a  satellite  about 
Venus  is 


P  =  1.6  hours 


(2.5.1) 


where 


a  -  638 7km 
/i*3.257xl0sfcm3/s2 
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The  mean  anomaly  is  a  fast  variable  with  a  short  period  (relative  to  the 
five  years  of  interest). 

The  argument  of  periapsis  <v  is  considered  a  long  period  term  and 
the  longitude  of  the  ascending  node  O  is  considered  to  be  secular 
(smoothly  varying  with  time)  term  (7:171). 

Since  the  fast  variable  M  is  just  an  oscillation  about  the  secular 
motion  of  the  orbital  elements  it  would  be  advantageous  to  average  out 
its  effects  on  the  overall  rate  of  change  of  the  orbital  elements.  This 
can  be  accomplished  by  setting 

(l  -  2 p  +  q)  =  0 


from  eqn(B.2c)  so 
<7=  2p-  l 

Satellite  Model 

The  satellite  model  is  concerned  with  three  parameters:  a)  the  mass 
of  the  satellite  m  (used  in  the  LPEf  the  solar  wind  model,  and  the 
atmospheric  drag  model),  b)  the  drag  coefficient  of  the  satellite  Ca  (used 
in  the  atmospheric  drag  model),  and  c)  the  projected  area  of  the 
satellite  A'  (used  in  the  solar  radiation  model  and  the  atmospheric  drag 
model). 

The  mass  of  the  satellite  consists  of  the  mass  of  several  satellite 
subsystems.  These  subsystems  are  the  Attitude  and  Control  System 
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(AOCS),  the  Telemetry,  Tracking,  and  Command  (TTfcC),  the  power 
system,  the  communication  subsystems,  the  satellite  antennas,  and  the 
satellite  structure.  Typical  values  of  these  subsystems  will  place  a 
nominal  value  for  the  satellite  mass  at  (20:67) 

m  -  1 085 Argr  (2.6.1) 

The  next  parameter  of  the  satellite  model  is  the  drag  coefficient  Cd- 
When  working  with  the  drag  coefficient  of  a  satellite,  several 
assumptions  must  be  made  concerning  atmospheric  molecules  (14:14-15): 

1)  The  satellite  is  considered  to  be  stationary  with  the  atmospheric 
molecules  flowing  past. 

2)  The  molecules  are  assumed  to  impinge  on  the  satellite,  be 
retained  temporarily  on  its  surface,  and  then  re-emitted. 

3)  The  collisions  between  incident  and  re-emitted  molecules  are 
neglected. 

There  are  several  factors  to  consider  in  calculating  C<r.  the  flow 
regime  through  which  the  satellite  moves  (reflected  in  the  Knudsen 
number),  the  mechanism  of  molecular  reflection  (presented  in  the 
accommodation  coefficient),  and  the  satellite’s  dynamics  and  orientation 
to  the  atmospheric  flow. 

The  type  of  flow  is  determined  by  the  Knudsen  number  Kn,  defined 
as  the  ratio  of  the  mean  free  path  of  atmospheric  molecules  to  the 
characteristic  linear  dimension  of  the  satellite  (11:184) 
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Kn  - 


X 

7 


(2.6.2) 


Assuming  a  Maxwell  distribution  of  particle  velocities,  the  mean  free 
path  is 

x-J~L 

where  o  is  the  collision  cross  section  of  a  molecule,  given  by 
a  •  4/rr 2 

n  is  the  molecular  number  given  by 

n  ■  — 
m 


with  a  simplified  density  model  for  the  atmosphere  given  by 


and  r  is  the  collision  radius  of  the  molecule.  The  Knudsen  number 
(eqn(2.6.2))  may  now  be  written  as 


Kn 


me 


« H 


/2(4  itr2)p0l 


(2.6.3) 


The  atmosphere  of  Venus  is  96.5%  carbon  dioxide  (6:173)  which  has 
a  collision  radius  of  (18:A.8) 


r-2.0xl0',2ftm 
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Other  typical  values  for  the  atmosphere  of  Venus  are  (18:157) 
m  ■  7 .3  1  x  1 0  ’29A:g 
a  *  4x  1 0  ~3km " 1 
H  «  1  10*m 

p0  »  3  .  !  9x  1 0'*kg/km 

and  assume  a  typical  characteristic  length  for  a  satellite  of 
/  -  5x  10  '3km 

These  values  yield  a  Knutsen  number  of  10.  In  order  to  have  free 
molecular  flow,  a  Knudsen  number  of  10  or  greater  is  needed  (11:184). 
This  implies  an  altitude  greater  than  110km  may  be  considered  free 
molecular  flow.  When  a  satellite  i6  in  this  flow  regime  its  drag  coefficient 
is  dependent  on  the  molecular  speed  ratio.  This  is  the  ratio  of  the 
satellite  speed  to  the  probable  molecular  speed. 

The  next  factor  to  consider  is  the  mechanism  of  molecular  reflection. 
The  energy  exchange  between  the  satellite  and  the  molecules  depends  on 
the  speed  and  direction  of  the  reflected  molecules.  It  is  assumed  that 
the  atmospheric  molecules  that  impinge  on  the  satellite's  surface  do  not 
reflect  specularly  but  attach  themselves  to  the  outer  layer  of  the 
surface  of  the  satellite  for  a  period  of  time  before  being  re-emitted. 
When  the  molecules  disassociate  from  the  satellite  they  are  emitted 
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diffusively,  having  lost  reference  to  their  original  direction  of  motion. 
This  diffuse  reflection  is  strongly  dependent  upon  the  nature  of  the 
satellite’s  surface  and  its  structure  (4:931). 

The  speed  of  the  re-emitted  molecules  is  a  function  of  their  kinetic 
temperature.  During  the  period  the  molecules  are  attached  to  the 
satellite  the  molecules  transfer  some  of  their  original  temperature  to  the 
satellite.  The  amount  of  uncertainty  (14:15)  in  the  new  temperature  of 
the  molecule  is  represented  by  the  accommodation  coefficient  o,  defined 
as 


a 


T,-Tr 

T,-Tt 


where 

Ti  =  original  molecular  temperature 

Tr  -  re-emitted  molecular  temperature 

Ts  =  satellite  surface  temperature 

King-Hele  suggests  (14:15)  that  the  accommodation  coefficient  is 
nearly  one;  therefore,  a  value  of  a  =  1  is  used  for  the  satellite  model. 
This  amounts  to  elastic  collisions. 

The  last  factor  to  consider  for  the  drag  coefficient  is  the  satellite’s 
dynamics  and  its  orientation  to  the  molecular  flow.  For  a  cylinder 
tumbling  end  over  end  the  drag  coefficient  can  be  related  to  all  the 
parameters,  above,  by  (4:940) 
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(2.6.4) 


1  =  cylinder  length 

d  =  cylinder  diameter 

With  a  =  1,  the  drag  coefficient  becomes 

2  (2.6.5) 

The  last  parameter  of  interest  for  the  satellite  model  is  the 
projected  area.  This  area  affects  the  atmospheric  drag  of  the  satellite 
and  the  pressure  of  solar  radiation  on  the  satellite.  The  size  of  the  area 
depends  on  the  orientation  of  the  satellite  as  it  orbits  Venus.  At  one 
extreme,  the  area  can  be  taken  to  be  the  projected  area  of  the  Bide  of  a 
cylinder  so  that  A-l+d.  At  the  other  extreme,  the  area  is  the  end  of  the 
cylinder,  or  A=xd2/\. 

It  would  be  unrealistic  to  assume  a  satellite  with  no  station  keeping 
abilities  would  remain  in  any  particular  orientation  with  respect  to  the 
atmosphere  or  the  solar  radiation.  This  precludes  the  use  of  either 
extreme  of  reference  area  mentioned  above.  An  uncontrolled  satellite, 
with  any  initial  rotational  motion,  would  begin  to  spin  about  its  axis  of 
greatest  moment  of  inertia.  This  is  due  to  relatively  small  external 
torques.  As  the  satellite  tumbles  about  its  orbit  it  can  assume  any 
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orientation  with  respect  to  the  planet's  atmosphere  and  the  solar 
radiation.  A  mean  value  of  the  projected  area  for  a  tumbling  satellite  is 
given  by  King-Hele  to  be 


A'  =  id(0.818  +  0.25d//)  (2.6.6) 

Typical  values  for  length  and  diameter  of  a  satellite  are 
1  =  11.60 

d  =  2.38 

These  values  give  a  projected  area  for  the  satellite  of 

A' m  24m 2  (2.6.7) 

To  summarize,  the  parameters  for  the  satellite  model  are 

m =  1085*g  (2.6.8a) 

C,-2  (2.6.8b) 

A '  ■  24m2  (2.6.8c) 
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III.  Approach  and  Results 

Several  simplifying  assumptions  were  made  in  order  to  model  the 
hypersurface.  The  first  assumption  was  with  respect  to  the  geopotential 
model.  The  other  assumptions  were  made  on  the  values  used  for  the 
orbital  elements. 

For  the  geopotential,  a  4x4  gravity  field  was  considered.  The 
coefficients  for  the  gravity  field  are  give  in  Table  3.2. 

Restrictions  on  the  orbital  parameters  consisted  of  two  types.  The 
first  type  of  restriction  was  to  choose  an  initial  value  for  a  parameter 
and  never  vary  it.  This  was  done  for  two  parameters,  the  mean  anomaly 
(M)  and  the  argument  of  periapsis  (o>).  The  mean  anomaly  is  a  fast 
variable  and  has  been  averaged  out  in  the  modeling  of  the  previous 
chapter.  For  this  reason,  any  initial  value  of  M  is  arbitrary.  Therefore, 
a  value  of  zero  has  been  chosen.  In  order  to  bring  the  hypersurface 
down  to  a  conceptually  understandable  four  dimensions,  another  constant 
orbital  parameter  was  needed.  The  argument  of  periapsis  was  chosen 
with  a  value  of  53'  so  calculations  of  the  time-to-impact  on  the  planet 
could  be  started  with  the  periapsis  of  the  satellite’s  orbit  in  eclipse 
whenever  the  longitude  of  the  ascending  node  of  the  satellite  is  zero. 
The  53*  relates  to  the  longitude  of  the  ascending  node  of  Venus  for  the 
epoch  used  in  this  study.  This  is  shown  in  Fig(3.1),  below. 
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Figure  3.1:  Orientation  of  the  Satellite  to  Venus  and  the  Sun 


The  second  type  of  restriction  was  in  the  range  of  the  parameters 
which  were  allowed  to  vary.  This  restriction  was  placed  on  two  of  the 
orbital  parameters,  the  eccentricity  (e)  and  the  inclination  (i).  The 
eccentricity  was  kept  small  (0.02<c  <0.10)  and  the  inclination  close  to  the 
equator  (0.5£i  315.5). 


i 
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This  leaves  two  orbital  parameters:  the  longitude  of  the  ascending 
node  and  the  semi-major  axis.  The  longitude  of  the  ascending  node  (fl) 
was  varied  throughout  its  entire  range  from  zero  to  2s.  The  semi-major 
axis  was  related  to  the  perigee  altitude  by 

hp  -  a(  1  -fl2)-  r  (3.1) 

and  treated  as  the  dependent  variable. 

Now  there  is  a  four  dimensional  surface  to  use  as  a  window  of 
survivability.  By  plotting  the  relationship  of 

/i,-F(i,e,fl)  (3.2) 

it  can  be  seen,  given  a  satellite’s  orbital  parameters,  how  low  in  altitude 
a  satellite  can  go  and  still  have  five  years  before  impact. 

Given  the  constraints  listed  above,  the  models  of  the  perturbations 
were  programmed  using  the  FORTRAN  programming  language 
incorporating  programs  from  JPL  (16).  The  computer  employed  was  an 
IBM  XT  compatible  computer,  specifically,  a  Connextion  XT  Turbo  with  an 
Intel  8088-2  microprocessor  and  an  Intel  8087  math  coprocessor. 

Once  the  program  for  integrating  the  LPE  was  constructed  a 
criterion  for  the  satellite’s  entry  into  the  atmosphere  was  needed.  In 
other  words,  how  low  in  altitude  must  the  satellite  be  before  it  is 

considered  to  be  within  the  atmosphere  of  Venus?  Recall  from  Section 

2.6,  the  interaction  of  the  satellite  with  the  atmosphere  of  Venus  can  be 
considered  as  random  molecular  collisions  for  altitudes  above  110km.  This 
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is  due  to  a  Knudsen  number  ^  10,  For  this  reason,  the  satellite  was 
considered  to  be  within  the  atmosphere  when  it  went  below  110km  in 
altitude.  A  sample  orbital  calculation  using  the  programmed  model  shows 
the  rapid  orbital  deterioration  when  the  periapsis  altitude  falls  below 
110km.  See  Fig ( 3.2 )•  below. 
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In  order  to  account  for  any  inaccuracies  in  the  model  a  margin  of 
safety  of  60  days  was  added.  The  satellite  was  considered  to  have  a  five 
year  survivability  if  it  went  below  a  110  km  altitude  after  5  years  and 
60  days  (1885  days). 

There  were  720  distinct  points  in  the  four  dimensional  space 
examined.  Each  computer  integration  run  of  the  LPE  took  approximately 
30  minutes.  It  took  roughly  four  or  five  runs  at  each  point  to  find  the 
initial  periapsis  altitude  needed  to  insure  a  five  year  survivability  of 
the  satellite. 

The  epoch  of  the  initial  conditions  was  23  February  1990.  The  data 
for  Venus  on  this  date  are  listed  in  the  table  3.1,  below  (2). 


1  Table  3.1:  Venus  Data  | 

Parameter 

Value 

Epoch 

19900223 

a(fcm) 

108208900. 

« 

.006777729 

«(°) 

2.6620179 

fl(°) 

53.555983 

a>C) 

-105.23265 

A/(°) 

43.540609 

n(Vsec) 

1.854317699e-5 

Prime  Meridian 

91.561209 
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The  coefficients  for  the  geopotential  are  listed  in  Table  3.2, 


|  Table  3.2:  Geopotential  Coefficients 

n 

m 

Cam 

Sum 

2 

0 

-.4520722821771336e-05 

. 00000000 00000000e+ 00 

2 

1 

•2039874602320317e-06 

.3981774131465276e+02 

2 

2 

.6582161633245868e-06 

-.61652959054 54911e+01 

3 

0 

.1342079593314583e-05 

.0000000000000000e+00 

3 

1 

.24962592008801 12e-05 

. 101474762506 3479e+02 

3 

2 

.269034 966 14638 19e-06 

.424676198132208e+02 

3 

3 

.3130468881843509e-07 

. 402251 9369940560e+02 

D 

a 

.24 134596 1 97689 15e-05 

. 0000000000000000e+00 

a 

D 

.703237661861 6599e-06 

.1321605322173613e+03 

4 

2 

.1151665471340579e-06 

.4338487639282226e+02 

4 

3 

. 1564751078966026e-07 

-.4872890600202502e+02 

a 

n 

.2868983692949129e-07 

•2042893824377979e+02 

Other  parameters  needed  for  orbital  calculations  are  given  in  Table 

3.3. 


Table  3.3:  Miscellaneous  Parameters 


Parameter 


Value 


#iC(*m3/8ec2) 


3.2485877e05 


Venus  Radius  (km) 

Venus  Rotation  Rate(' 
/sec) 

*it„C(*mVsec2) 


6051.0 

-1.71460706e-05 


0.13271244el2 
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The  next  step  in  the  development  of  the  hypersurface  was  to  fit 
the  data  generated  from  the  model  to  a  least  squares  curve.  The  first 
step  in  deriving  a  least  square  representation  in  multidimensions  was  to 
define  the  dependent  and  independent  variables 

Xi  =  independent  variable(s)  (may  be  multidimensional) 

yi  =  measured  dependent  variable 

For  a  least  squares  approximation  the  data  is  to  be  fitted  to  a 
multiple  of  functions 

y  •  c,/,(x)+  c2/2(x)+  ...♦  cn/m(x) 


or 


(3.3) 

i-i 

Assume  that  there  are  n  measurements  of  y,  so  eqn(3.3)  may  be 
expanded  for  each  value  of  y. 

yt-  c1/1(x1J+c2/2(x,)...c*/m(xI) 

y2-c1/i(x2)*c2/2(x2)...c,n/in(x2) 

y.“C,/,(x.)*c2/2(x„)...cw/m(xll) 
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Define  the  residuals  as 


m 

*im  y  i  - 


i-l 


m 

^2-y2-Lc,/,(x2) 


cl 


m 

/?„-y Xc«/«(x») 


<*i 


In  order  to  get  the  best  least-squares-fit  to  the  chosen  functions,  the 
sum  of  the  squares  of  the  residuals  mu6t  be  minimized.  That  is  to  say 


n 

R~Y.R2im  minimum 


(3.4) 


/-I 


Define  the  total  derivative  of  the  residual  (for  later  use)  as 


dR  d*  dR 

dR  •  - — ac,  ♦  - — dc,+  ...♦  - — dc_ 

ac,  ac j  ac* 


To  find  the  minimum,  set  the  derivative  of  eqn(3.4)  to  zero. 


**.0. -Lf  R>m2f  R  — t 

*ct  act/4i  1  ;ac» 


*  -  1,2,3 . m 
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dR.  ,  s 

but  jt:"'*1*') 


80 


n 


I 


*,/.(*,) 


Z  y/'Zc^'(xi)  /.(x<)-0 

/-  L  <-> 


(3.5) 


Rearranging  terms  in  eqn(3.5)  yields 


£y ;/*(*;)“  Zc<  £/<(*/)/»(*/) 

/-i  i-i  /-i 


*-1,2 . m  (3.6) 


Next,  put  eqn(3.6)  in  matrix  form  for  ease  of  manipulation.  In  order  to 
accomplish  this,  some  definitions  need  to  be  made.  Define 


/-i 


«*•■£//(*;)/*(*/) 


so 


c  ,a  kl 


Now,  define 
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/7i(*i)  /2(*J 

I  /i(*a)  /2(x2) 


/ Ax J 


so 

j4  -  F  rF 
B-  Fry 

Therefore,  eqn(3.7)  may  now  be  written  as 

Fry«(FrF)C  (3.7) 

Since  y  is  given  and  F  is  chosen  to  fit  the  data,  the  interest  is  in 
finding  a  C  that  will  minimize  the  residuals.  This  will  give  a  best  fit  of 
the  data  for  the  chosen  functions.  The  vector  C  is  (from  eqn(3.7)) 

C-(FrF)~'(Fry)  (3.8) 

It  now  becomes  a  matter  of  judicious  selections  of  functions  to  best 
fit  the  calculated  periapsis  altitudes. 

The  calculated  values  for  periapsis  altitudes  for  inclinations  of  0.5' 
and  1.5*  are  shown  in  Table  3.4  (at  the  end  of  this  section)  along  with 
graphs  of  the  results,  Fig(3.1)  and  Fig(3.3),  and  contour  maps  of  the 
data,  Fig(3.2)  and  Fig(3.4).  The  remainder  of  the  calculated  data  may  be 
found  in  Appendix  D  along  with  graphs  in  Appendix  E. 
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From  analysis  of  the  data  the  periapsis  altitude  (hp)  behaves  as  a 
sinusoidal  function  with  the  ascending  node.  The  hp  behaves  as  a 
inverse  fourth  root  in  eccentricity  and  as  a  quadratic  in  inclination. 
See  Fig(3.1)  and  Fig(3.3).  For  the  trends  in  the  longitude  of  the 
ascending  node  refer  to  Fig(E.33).  The  functions  used  to  fit  the 
calculated  periapsis  altitude  are 

/,({.«. n)-e‘,/4 

sin(rt-  149*) 

/,(!.•. fl)-i(l  -i)/r/180# 

These  functions  result  in  a  least  squares  fit  equation  of 

hp-  1 00 .6829a*l/*  +  48 .7606sin  (D  “  1 49°  )♦  0. 10607  i(i+  103.9) 

where  the  inclination  is  in  radians.  The  correlation  coefficient  for 
eqn(3.9)  is 

r  =  0. 97257  (3.10) 

where  the  correlation  coefficient  is  defined  by  (10:87) 


with 
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a 


kp  .x 


Y.ih*-hp«y 

e-i 

nT 2 


f  I 

/l-  1 

'*  i- 1 

hpi  =  calculated  value  of  hp  from  the  model 
hpx:  =  least  squares  calculated  value  of  hp 

The  least  squares  representation  stays  within  ±20km  throughout 
most  of  the  range  of  parameters  examined.  This  is  a  10%  variation,  well 
within  engineering  practice.  The  range  of  data  that  does  not  support 
this  small  variation  is 

10.5°  <i <  15.5° 

e=  0.02 

-30°  <n<60° 

In  this  range,  the  error  is  as  high  as  50km,  however,  all  the  deviations 
from  the  computed  model  are  on  the  positive  side  of  the  hypersurface 
and,  therefore,  conservative.. 
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This  region  of  high  deviation  constitutes  less  than  3X  of  the  entire 
hypersurface.  If  it  is  excluded  from  consideration,  then  the  least 
squares  equation  becomes 

hp -  100.73e'1/4+  50.265sin(fl-  1 49#)+  0 .0872j(i  +  1) 
with  a  correlation  coefficient  of 

r  =  0.98570 
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Flo  3.5:  Per  lops  Is  Altitude  (  l  =  1.5  deg) 


Fig  3.6:  Periapsis  Contour  (i  =  1.5  deg) 


I tude  (  l  =  0.5  deg) 


Fig  E.8:  Periapsis  Contour  (i  =  0.5  deg) 


Conclusions 


IV.  Conclusions  and  Recommendations 


This  effort  investigated  the  creation  of  an  analytic  function  in  the 
six  dimensional  orbital  parameters  space  that  separated  a  satellite’s  orbit 
between  a  five  year  or  greater  survivability  and  a  satellite’s  potential 
for  entering  a  planet’s  atmosphere  in  under  five  years.  On  the  basis  of 
this  investigation,  the  following  conclusions  are  made: 

1.  It  is  possible  to  generate  smooth  curves  in  the  orbital 
parameter  space  that  separate  a  satellite’s  life  expectancy  as  a 
function  of  periapsis  altitude  and  time. 

2.  The  periapsis  altitude  ( hp )  behaves  as  a  sinusoidal  function 
with  the  ascending  node.  The  hp  behaves  as  a  inverse  fourth 
root  in  eccentricity  and  a  quadratic  in  inclination. 

3.  It  is  possible  to  generate  a  survivability  function  hp  for  the 
hypersurface  with  an  acceptable  engineering  error  of  10%.  The 
approximation  generated  in  this  thesis  is 

h„-  1 00 .6829  e’ 1/4  +  48 .7606  sin[(/2  -  1  49#)rf/ 1 80]+  0. 10607i(i  +  1) 
with  a  correlation  coefficient  of 

0-0.98570 
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where  the  constraints  in  thi6  thesis  are 


M- 0° 

co  -  53° 

0.02  <e  <0. 10 

0°  <i<  15° 

0°  </2  <360° 


Recommendations 

Based  on  the  findings  of  this  investigation,  the  following 
recommendations  for  further  study  are  proposed: 

1.  Investigation  of  the  behavior  of  the  satellite  for  inclinations 
up  to  a  value  of  90"  and  eccentricities  above  .10  should  be 
conducted. 

2.  Further  study  of  the  modeling  of  the  hypersurface  using 

simpie  functions  (sin  x,  cos  x,  ex,  etc.) 

3.  A  layered  mapping  of  five  year  increments  above  the  planet 
Venus  should  be  performed  so,  should  a  satellite  lose  its 
station  keeping  abilities,  a  determination  can  be  made  at  once 
concerning  the  lifetime  of  the  satellite. 
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Appendix  A 


Conversion  of  the  (W>pr>tpntriAl  into 
the  Classical  Orhital  Elements 

Although  the  following  derivation  has  been  outlined  by  Kaula 
(13:30-37)  it  is  included  here  for  clarification  of  the  geopotential  model 
in  Section  2.2.  This  derivation  of  the  geopotential  deviates  from  that 
presented  by  Kaula  in  that  Hansen's  coefficients  are  used  in  formulating 
the  Eccentricity  Function. 

F(r,0,*)- />"(sin©)(C ,mcosmj  +  S linsin/n^)  (A.l) 

The  geopotential  in  its  spherical  harmonic  representation 
(eqn(2.3.15))  must  be  translated  into  a  Keplarian  representation.  To  do 
this,  two  special  functions  will  be  introduced:  the  Inclination  Function 
( Fnmp{ i) )  and  the  Eccentricity  Function  {Gnpg(e))t  both  to  be  defined 
later.  In  order  to  use  these  two  functions  the  geopotential  must  first  be 
rearranged. 

The  first  step  will  be  to  get  the  geopotential  into  a  form  compatable 
with  the  Inclination  Function.  The  following  relationships  will  be  needed. 

cosmx-  Ro[ Rq { ( cos x ♦  ysinx)*1} 
or 


66 


cos  mx 


-*«{i(r)''cos"‘xei"'x} 


In  like  manner,  note 


sin  mx  -  Ro 


{Um*)r'c°*'"xsiD'x} 


Multiply  the  two  equations,  above,  to  get 


sm  jccos 


,  ( o,x  -  Q'txY ( a1**  V 

X"l  2y  j  (  2  j 


*f>/*  Q~*iX 


My 


-K)l*Q-klx 


sin 


x{cos(a  +  0  - 2g - 2h)x  *  j  sin(a  +  0 -2g-2h)x} 


also 


cos  a  cos  0  -  -[co8(a+  0 )♦  cos(a-  0)3 

m* 


{A. 2a) 


(/1.2b) 


(/4.2c) 


(/1.2i  ) 
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sin  a  sin  0  -  -[- cos(a  ♦  0)  +  cos(a-  0)] 


(.A  2a) 


sin  a  cos 0  •  -[sin(a  +  0)+  sin(a-  0)] 


(A.2f) 


cos  a  sin  0m  -[sin(a+  0)-  sin(a-0)] 


(A.2  g) 


cos  (a  *  b)  =  cos  a  cos  b-  sin  a  sin  b 


(A.2h) 


sin  (a  ♦  b)*  sin  a  cos  b*  cos  a  sin  b 


(A. 2i) 


And,  some  spherical  trigonometry  identities  that  we  will  use  are 


cosa  *  cosbcosc*  sin bsin ccos  A 


(.A.2  j) 


cosb  *  cos  ccos  a  +  sin  c  sin  a  cos  B 


( A.2k ) 


cosc*  cos  a  cos  b+  sinasinbcosC 


(A. 21) 


sin  A  sin  B  sin  C 
sina  sinb  sine 


(A  .2m) 


Now  the  geopotential  is  in  a  form  compatible  with  the  Inclination 
Function.  The  geopotential  may  be  written  as 


K-n>, 


(A.3a) 


where 
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(A. 3b) 


uR1 

V  imm  — (sin  0)(C  lmcosmf  +  S,„sinmt) 
r 

Let 

m^  -  m(a  -  0#)  +  mCl-  mCl 

where  Qtm  prime  meridian  of  Venus 

a  =  right  ascension 

Cl  =  longitude  of  ascending  node 

or 

mQ  m  m(<i~  fl)+  m[fi-  Q t) 

Using  eqn(A.2i)  and  eqn  (A.2j),  we  have 

cos mj>  •  cos[/n(a  -  fl)]cos[m(/2  -  ©J] 

-8in[m(a-  fl)]sin[m(fl  -  0#)] 
sinm^-  sin[m(a-  fl)]cos[m(fl-  ©,)] 

+  cos[m(a-  fl)]sin[m(fl-  6,)] 

Using  eqn(A.2k)  and  eqn(A.21),  we  get 
cos(co  ♦  /)  -  cos(a  ♦  /i)cos© 

cos©  *  cos(cu  ♦  /)  cos(a-  Cl )*  sin  (to  +  /)sin(a-  /3)cosi 
Eqn(A.2m)  yields 


(A. 4) 

(A. 5a) 

(A. 5b) 

(A. 6a) 
(A. 6b) 
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sin  0  =  sin(a>  +  /)sin  i 


(A. 7) 


Eliminating  the  terms  containing  cos(to+f)  from  eqn(A.6a)  and 
eqn(A.6b)  gives 


cos0  -  cos2(a  -  /2)cos0  +  sin  (to  ♦  /)cosisin(a  -  fl) 


sin(a  -  /})  -  sin(<o  ♦  /) 


cost 
cos  0 


(A.  8) 


Using  eqn(A.2a)  and  eqn(A.2b)  for  (a  +  fi)  in  eqn(A.5a)  and  eqn(A.5b) 


gives 


cosm 


(  m  /  _  \ 

0“/?e<^(  )/‘cos,n't(a  - /2)sin '(a  - /2) 

(t-oV  s  J 


x  [cos(m(/2- 0 J)+ ysin(m(/3- 0r))]}  (A. 9a) 


sin  mfl  -  /?o|  ^  jy*c°sm**(a  -  /3)sin’(a-/3) 

x  [sin(m(/3-0#))-ycos(m(/2-0j)]}  (A. 9b) 

Applying  eqn(A.6a)  and  eqn(A.8)  to  the  above  equations,  eqn.s(A.9), 
yields 


cosm^ 

sinml 


/)  8in*(to+  /)cos'» 


cos"© 


[cos[m(/l-0  ,)]♦  >sin[m(/2  -  eJjl'N 

8in[m(n-0j]-yco8[m(/2>eJ]JJ 
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(^.10) 


Now,  substituting  eqn(A.7)  and  eqn(A.lO)  into  eqn(A.3),  using 
eqn(A.12c),  results  in 


Jm 


nR‘,ir 
•777 ;Lr 


tint 


sin,'m'2t(0 


c  0 


xRo[[Ctm- jSlm]cos[rn({l-et)]+[Slm+  /Clm]sin[/n(fl- 0,)] 

x  X  f  m  j/‘cos'r,'‘(to  -  /)sin‘‘"’2t*‘(£o  +  /)cos‘(t)l  (>4.11) 


where  A:  =  integer  part  of  (i-m)/ 2 
In  eqn(A.2c),  let  a=l-m-2t+s  and  fi=m-s,  then  define 

4  •  sin  l  n  2,'c((v  +  /)cosm  s(oj  +  /) 

2 1  *  /r0  V  0  J\hJ 

x{cos[(l-2/-2g-2A)(oo  +  /  )]«•  /sin[(f-  21-  2ff-  2fc)(o>  +  /)]>  (>4.12) 

Insert  eqn(A.12)  into  eqn(A.ll): 


where 

Am\C tm  -  jS llB]co8[m(fl- 0#)]+  [5<(n+  jC lm]sin[m(/3-  0#)](>4. 13b) 
For  notational  convenience,  let 
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6  ■  m(/2  -  fl#) 


k*(l-2t-2g-2h)(u)  +  /) 


(A. 14a) 
(A.l4b) 


Now,  using  eqn(A.2d-g),  we  have 


A{cos[(l-2t-2g-2h)(a>+  /)]  +  jsin[(l- 2t-  2g-  2h)(co  *  /)]> 


-  [C/m  -  /5,m][cos6cos k*  j cos 6 sin  A.] 
*f-s’/in~yCz».][siii6cosA.+  7sin6sinX] 

[Clm~  jS  ,m]cos(6+  k)+[Slm+  jC  lin}sin(6+ k)  {A.  15) 


Inserting  eqn(A.15)  into  eqn(A.13b) 


IS  (m  iVi  /  .T  Imt  SiD  (0 

r  r*o 


x*.  j.nc-io^-cr- . ‘v . 

,  •  o  V  s  /  2 


x . rrv'*sxv>->>' 


x  [clm- jS,m]cos{6+  *.)♦[$,„*  ;Ctin]sin(0+  A))  (A. 16) 


Now  we  will  eliminate  the  complex  numbers  in  the  summations.  Note 


If  (ha)  is  odd  and  k  is  the  integer  part  of  (i-m)/ 2,  then 


-  (-  1  -/(-  1  )"*--/(-  1  )*“ 


(A. 17b) 
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However,  if  ( l-m )  is  even,  then 


j‘(- -  (-  1  )*’*-(-  1)'** 
Combining  eqn(A17b)  with  eqn(A.17c)  yields 


f  -  1  )'**  ,  l-m  even 

\(-  1)*'*  ,  l-m odd 


Eqn(A17d)  shows  the  real  part  of  the  geopotential 
depend  on  whether  (l-m)  is  even  or  odd.  In  other  words: 


(mf Sin  (I) 


f-0 


(Clmcos(6+  A)  +  5\msin(6+  A)  ,  l-m  even 
-5‘lincos(6+ A)+ C)msin(6+ A)  ,  l-m  odd 


where  6+A-(f-2f-20-2/i)(cu+/)+m(/2-*0j 


Let  p*t+g+h  and  define 


H  ■  6+  A  -  (/  -  2p)(w*  f)+m[{l-Qt) 


Using  eqn(A18b) 


i  /  ^  ^  (  \  »  f  |  v  fc*  I  i 

v imm  -77rLT««t(“  1)  Sin  (0 
r  ,.o 


(i4 . 17c) 

(d. 17d) 

(eqn(A16))  to 

-  1)' 

(if.  18a) 

(A.  186) 
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X 


my&n, . rt*r-')(r)<-'>- 


*(C ,mcosH  *  S lmsinH  ,  l-m  even 
\- S ,ncosH  +  C,„8in H  ,  l-m  odd 

The  summations  of  the  geopotential  (eqn(A.19))  must  now  be 
interchanged.  To  do  this  we  must  note  some  relationships.  If  t< 0  then 
Ti»*=0  (since  I  (-*)  I-  °°).  Also,  if  t>k  then  Th&O  (since  |(  i- m  -  2t)  |- *° ). 
Hence,  the  sum  on  t  can  be  -«<!<•,  since  Ttn  will  be  zero  for  all  values 
of  t  except  0<t<k.  Also,  the  summation  on  p  can  be  -®  <  p  <  ®,  since  the 
binomial  coefficient 


m  -  s 

p-t-g 


) 


will  be  zero  if  p<t+g  or  p>t+g+m-s.  Therefore,  the  limits  on  the 
summations  are  independent  of  each  other  and  may  be  interchanged.  So, 
substituting  for  Timt,  eqn(A.12c),  the  geopotential  (eqn(A.19))  may  be 
written  as 

y  ,.***!  f  y  (2(-20!sin(,""'2,)(Q 
r'*1  !(!-!)! (I- m-2«)f2a,’a* 


(C  ,„cosH  *  S  Jm8in  H  ,  l-m  even 
-  S  ,mcosH  *  CtmsinH  ,  l-m  odd 


(A  .20) 
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Before  the  Inclination  Function  can  be  introduced  into  the 


geopotential  (eqn( A20) )  the  limits  of  the  summations  over  p,  t,  and  g 
must  be  evaluated.  The  minimum  of  p=t+g+h  will  be  zero  since  the 
minimum  values  of  t,  g,  and  h  are  zero.  The  maximum  value  of  p  will  be 

(P)»a«"(-*+  0.„- l  (/.21a) 

since  1 2.  0.  Also, 

(0 .«-(p-0-*).„-  P 

or,  if  p>k  then  as  can  be  seen  from  eqn(2.2.12c)  for  7W 

Therefore 


(O 


max 


Pik 
p>  k 


(A.  21b) 


Now  to  look  at  the  index  g  in  eqn(A.20).  The  range  of  g  will  be 
limited  by  the  binomial  coefficients  involving  g.  So,  look  at 

( l-  m-2t+  s  V  m-  s  \ 

l  g  Jlp-f-eJ 


The  above  binomial  coefficients  have  the  following  four  conditions: 

(9*  0 

m-  st  p- 1- g 
l-  m-  21*  s  £  g 
p-tig 
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The  first  two  conditions  give  the  lower  limit  for  g  and  the  last  two  give 
the  upper  limit  of  g.  These  limits  are 


gi«max[0,p-<-m+s] 


The  Inclination  Function,  Fiapli)  (7:110),  is  defined  as 


(  2  /  -  2  / )  ■ 


_  ...  Vs  (<£1-41 


2<)l  2 


2l-2t 


sin  c,-"-ao(«) 


(yl.21c) 

( d  .2 1  d) 


(A. 22) 


Inserting  the  Fimpd )  (eqn(A22))  into  the  geopotential  (eqn(A20)) 
yeilds 


V 


r1*1 

r  ,-o 


Imp 


(0 


{C(incos// +  5“(insin //  ,  I-m  even 
-5‘(ncos//  +  C,msin//  ,  l-m  odd 


(A.23) 


The  geopotential  (eqn( A23),  is  now  ready  to  be  modified  to  allow 
the  inclusion  of  the  Eccentricity  Function,  GWe).  The  portion  of  the 
geopotential  to  be  transformed  is 


1  (C lmcosH  +  S ,msinH  ,  l-m  even 
S  lncosH  *  C  lmsin  H  ,  l-m  odd 

To  do  this,  examine 
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(A.  24) 


1  (cos H 
r^'lsin  H 


Only  long  period  terms  (those  which  do  not  contain  M)  are  of 
interest.  So,  the  geopotential  (eqn(A.23))  may  be  averaged  with  respect 
to  M.  To  do  this,  integrate  with  respect  to  M  from  0  to  2*  and  then 
divide  by  2*.  This  is  the  Method  of  Averaging  (17:11).  The  following 
identities  will  be  used 

...  r2df 

d  M-—.*—  (A. 25a) 

a^l-o2 


r.  ai. 

1  +  acos / 


(A. 256) 


r 


a(  1  -  a2)  a(  1  -  a2) 


cos  / 


(A. 25c) 


Prom  eqn(A.25b)  and  eqn(A.25c)  deduce 


The  binomial  expansion  of  (l+ecos /)w  is 

( 1  +  acos/)1'1  ■  ^  *  je'^cos1*1/ 

Integrating  eqn(A.24): 


(A.25d) 


(A. 25/) 
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J_  f  2‘  1  /  cos  H 
2/rJo  r^'lsin// 


dU- 


1  ram  r2  (cos  H 

2/rJo  r'*'a2(  1  —  ©a)'/2\sin// 


2  n 


C- 

■'o  ( a 


(  1  ♦  acos /) 


l-\ 


(aJVl  -  o 2 )  [  a  ( 1  -«2)] 


/costf 
‘-'Isin  H 


1  f2*  (l  +  acos/)*'1  (cosH 
2/rJo  a'*'(l -02)''1/2\sin//  d/ 


df 


cos  H 
sin  // 


d/  (>1.26) 


Here  the  long  period  terms  are  of  interest.  This  is  when  (see 
Section  2.5) 


(J-2p)±(fe-2d)«0 

Recall  H  is  given  by  eqn(A.18b).  Expand  the  sin  H  and  cos  H  using  the 
identities  of  eqn.s(A2/i)  and  (A2j)  and  let 

fi-  (1-  2p)co  -  m(/l-  0) 


so 


sin[(/-  2p)/+2?]«sin(/  -  2p)/cos8  +  cos  (/  -  2p)/sin  B 
cos[(l-  2p)/  +  B]»  cos (i-  2p)/ cosB -  sin(i-  2p)/sin  B 
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The  sine  function  is  odd,  and,  therefore  (over  the  limits  of  integration) 
may  be  ignored.  So 


sin[(/~  2p)/  +  2?]«  cos (i  -  2p)/sin  B 

cos[(Z-  2 p)f*  S]=*  cos (i  -  2 p)/ cos B 

Inserting  these  identities  into  eqn(A26),  then  pulling  out  variables  not 
dependent  on  /,  and  interchanging  the  integral  and  summation  results  in 


— 1 _ vf'-'V 


x 


X2^Z  cos4/ cos  (  l-  2  p)f  d  f 


sin(Z-  2p)a>+  m{Q-  6) 
cos(/-2p)a>  +  m ( H -  0) 


(yi.27) 


In  order  to  integrate  eqn( A27)  Hansen’s  Coefficient  will  be  introduced. 
Hansen’s  coefficient  may  be  defined  as  (2) 

3/l(-l)b(n  b+1)*k2^r/0  cosb/cos mfdf  (yl.28) 
let 

n« -{l*  l) 
m  -  (l  -  2 p) 

Then 


*(!♦ 
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using(l:ll) 


') 


So 


v-((- l).(l-2p) 
A  o 


(l-«2) 


cos h  fcos(  l-  2  p)/d/ 


Insert  this  Hansen’s  Coefficient  into  eqn(A27)  and  the  integration 
becomes 


J_  (2‘(sinH 
2rtJo  \cos  H 


dM  - 


1  y.(<«D.«-a»)/8in[(l-2p)/  +  m(n-g)] 

"a'*'  0  \cos[(/-2p)/  +  m(/2-0)] 

or 

_j_  /sin[(i- 2p)/  + 

“  a  *^(2^-0 \cos[(i-2p)/  +  m(/2-0)] 

Where 


r  m  y-(l*l>. <»-**> 

uIK2p-D  "  *  0 


28a) 


28b) 


I 

I 


to  follow  the  notation  of  Kaula  (12:110). 


Comparing  eqn(A.28)  and  eqn(A.28a)  shows  n=-(l-l).  Recall  from  eqn(A.l) 
that  1  £0;  therefore,  n<2  and  Hansen's  Coefficient  for  n^2  may  be  UBed 
which  is  defined  as  (2) 


y  ft 
^  0 


flVm'  1  y  Y  ~n~2  Vl'ni  +  2dV®V 

2  J  h>  [lml  +  2dj[  d  ){2J 


let 


n- -(Z+  1) 


m  -  /  -  2  p  ^  |  m 


'  p'm  p  ,  p£  1/2 
>p'm  l-  p  ,  p>  1/2 


(  ^  .29) 


So,  Hansen’s  Coefficient  becomes 


> _ ! _ ?(  )(2d+l- *>')(*-) 

(l  -  q2)1’''2  &o\2d+  l-2p')\  d  JK2  ) 


and  the  Eccentricity  Function,  GipOp- u(e),  is 


• . «-(r*-SU.VX“T''») .  <-■“> 


So  that  the  integration  of  eqn(A.26)  is 


jl  r2*j_/8in//  dW _ 1 

2  ir  Jo  r  \co8// 
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'r-7  l'l  V2d+  1  -2p')(  oY‘"~2p  fcosH 
£0Ud+ i-2p'](  d  2  J  (sintf  (  '  } 


Note  that  the  factor  of  1/2  in  the  summation  over  d  has  been  dropped. 
This  is  because  the  two  terms  satisfying  the  long  period  variation, 
(l-2p)*(b-2d)- 0,  are  symmetric  in  the  binomial  expansion  and  combine  to 
make  the  (2 d+l-2p)  substitution  for  b. 

Now  the  disturbing  function  for  a  nonspherical  planet  is  given  by 

»  ) 

{A. 32a) 

1-2  tn-0 


i  i 


where  V  ,m  -  £  F  lmp(i)  £  C,M(o)S<in,#(a> ,  hi  .  fl . 0)  (yl.32b) 

a  pT0  qVTm 


a 


/Cllncos$+ S,m8in*  ,  i-n 
n  \-S,„cos^+Cfmsin^  ,  £- 


-  m  even 
m  odd 


( A .32c) 


where  0  -  (l- 2p)u>  *  (l  -  2p  +  q)M  *  -  0) 


(A -32d) 
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Appendix  B 

Conversion  of  the  Geo  potential  into 
the  Modified  Qrhit.nl  Elements 


In  Appendix  A  the  geopotential  was  found  to  be 


-  i 


1-2  m-0 


where  V  lm 


hr: 


7rfZf^(0  I  ClM(05lBM(fl),tf,fl,9) 


p-o 


and  S lnf>t 


1C ,ncosj+  Slmsint  ,  l-m  even 
\-S lncos1>+  Cl(Bsin0  ,  l-  m  odd 


where  $  =  (l- 2  p)a>  *  (l- 2p  *  q)hi  +  m(Cl- Q) 
Define 


then  V*  (eqn(A326))  may  be  written  as 

a  p*o  ♦— 


(A .32a) 

(A. 326) 

(A. 32c) 

(  A  .  32<J ) 


(B.la) 
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where 


$  Impt  ( to  ,  A/  ,  Cl ,  9 ,  k  N 


(cos^‘  ,l-m  even 
sin odd 


(B.lb) 


t‘-(t-2p)co*(l~2p+q)A/  +  m(n-0)-mXlm  (B.lc) 

For  Venus  only  the  long  term  effects  are  of  interest  so 


q  =  2p  -  l 


as  stated  in  Section  2.6.  Replace  S**pq  with 

.  _  /cos(f - quj).  1- m  even") 

“"f  \sin(f - <7<u), l- m  odd  J 

( B.2a ) 

where 

(8.2t>) 

q*  2p- l 

(B  .2c) 

m  =  0,1,2 . 1 

(B.2d) 

Multiply  and  divide  eqn(A2)  by 


e 


1*1 


and  use  the  trigonometric  identities  of  eqn.s(A2h)  and  (A.2i),  to  yeild 


S 


Imf 


_l_  fcosf  e'^cosqcu  ♦  sinf  e'^sinqco ,  l-m  even 
•  ,f,V8inf  e"lco8<7o>  -  cosf  elflsin  <jto  ,  l-  m  odd 


(*3) 
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and  define 


(«•  4) 


Note  f  (eqn(B.2b))  is  a  function  of  X*.  Also,  note  that  •l',«inq«>  and 
are  related  to  h  and  k.  From  eqn(A.26) 


,flsinqco  •  j(>)‘  '  cos,"‘o>sin'tt>  j 


But  can  be  written  as 


0m.  a'«i‘V 


so 


'•'sinqco  -  j  (/)''« l,1''cos*',co®‘sin'coJ 


Recall 

h  *=  esin  co 

k  ~  0  cos co 


BO 


©'’'sin 


(2.1 .2a) 
(2.1 .26) 


To  extract  the  real  part  of  the  above  eqn,  let  i  =  2>3  so 
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J*  I 


Q  V-/  1  \l(  Ifll  li.l«l-2i-3|,2i-3 


(5.5a) 


where  n, - |q  |  +  3)/2 


Similarly,  •l*,co*qa>  may  be  expressed  as 


C0SQ60-  2^(-  1  )'^  .  Jfc  « 


2; 


(5.56) 


where  n2  =  \q\/2.  So  that  the  geopotential  becomes 


Vln-^J,nY.F‘> «p(‘)  I  G',M(s)5Imf(h,*.X„)  (5.6a) 


<  i 


p-0 


where 


c  »ol,l5‘’ 

J  Imt  9  J  lm» 


(5.66) 
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When  calculating  Lagrange’s  Equations  for  the  Geopotential  in 


Section  II  the  following  derivatives  are  useful(lla). 


dS 


Imq 


dh 


dS 


dk 


dS 


Imt 


ts I  q  dS2  .  . 

cosf  +  m  even 


dh 


I  Q  i  4ft 


dSt  q  dS2 

-sinf  +  rr7TT"c08f  •  m  odd.. 


dh 


q  I  4ft 


dS 2  q  dS i  . 

— r-cosf  +  r-^-r-r —  smf ,  l- m  even 
dh  |  q  I  dh 


dS  2  q  dS  | 

——sin  f  ♦  ~cosf ,  l  -  m  odd 

dh  <7  4ft 


-  m < 


-5, sin f ♦  — ^-Sjcosf  ,  l-m  odd 
I  <7  I 


dQ.  1  r.  „  <7  -  .  „  , 

$,cosf+- — -5,8inf  ,  i-/n  even 


(C.la) 


(C.  lb) 


(C.lc) 


lac'  2(l-l/2)„^  _1 _ 

;17"  l-eJ  C  2,-2>'*,!(1-82),-,/j 


(C  .2a) 


(C  .2b) 


XV' _ 

t4l(2p'-2fc 


fc(*-l)? _ fiV 

-  I  )!(fc  ♦  I-  2p')\  fc  !  v 2 > 


(C.3) 


f  (2i-2t)»(/-m-20 

- -«COSl) - - - rr^;8in  I 

di  ht\(l-t)\{l-m-2t)\2il  2t 


WrK-'K 


(2/-2Q! _  . 

8in,,4o<!  (i  -  0-  2f)^  22I*2,SI 


.  (J-m-Jt). 

sin  i 
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Aprendis  £ 


Supplemental  Tables  far  Sectfon  HI  And  IY 

On  the  following  pages  are  tabulations  of  the  data  for  inclinations  from 
2.5  degrees  to  15.5  degrees.  The  data  ranges  are: 

2.S0  <i<  IS. 5° 


0.02  <e<0. 10 

0°  </l<330° 

co -53° 


A/  -  0 


i 
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TABLE  D.l:  INSERTION  WINDOW  HYPERSURFACE  (<-2.5*  .  3.5*) 


«(*) 


to 

e 

no 

Hp  (Km) 

2.50 

0.02 

0 

2.50 

0.02 

15 

2.50 

0.02 

60 

2.50 

0.02 

105 

2.50 

0.02 

150 

2.50 

0.02 

195 

2.50 

0.02 

240 

2.50 

0.02 

285 

2.50 

0.02 

330 

2.50 

0.04 

0 

2.50 

0.04 

15 

2.50 

0.04 

60 

2.50 

0.04 

105 

2.50 

0.04 

150 

2.50 

0.04 

195 

2.50 

0.04 

240 

2.50 

0.04 

285 

2.50 

0.04 

330 

2.50 

0.06 

0 

2.50 

0.06 

15 

2.50 

0.06 

60 

2.50 

0.06 

105 

2.50 

0.06 

150 

2.50 

0.06 

195 

2.50 

0.06 

240 

2.50 

0.06 

285 

2.50 

wmm 

330 

2.50 

0.08 

0 

2.50 

0.08 

15 

2.50 

0.08 

60 

2.50 

0.08 

105 

2.50 

0.08 

150 

2.50 

0.08 

195 

2.50 

240 

2.50 

0.08 

285 

2.50 

0.08 

330 

HEU!!] 

0.10 

0 

2.50 

0.10 

15 

2.50 

0.10 

60 

0.10 

105 

2.50 

0.10 

150 

0.10 

195 

Ikl i* 

0.10 

2.50 

0.10 

285 

2.50 

0.10 

330 

241.48 
229.82 
209.24 
228.35 
265.59 
295.09 

305.57 
294.30 

265.49 


195.14 
182.38 
165.19 

190.92 

236.14 
270.89 

280.58 
262.34 
223.94 


173.68 

161.74 

146.98 

172.46 

217.95 

253.49 

262.23 

241.92 

202.07 


160.56 

149.43 
136.09 
160.38 
204.54 

239.41 

247.41 
226.53 

187.43 


151.17 

140.82 

128.58 

151.53 

193.92 

227.49 

234.69 

213.99 

176.55 
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TABLE  D.3:  INSERTION  WINDOW  HYPERSURFACE  (i-6.5*  .  7.5*) 


«•)  e  no  HP  (Km)  I  *(•)  e  fl(°)  Hp  (Km) 


260.98 

250.40 

223.55 

227.17 
260.10 
292.93 
308.61 
303.22 

281.17 


202.73 

190.54 

170.76 

189.00 

231.14 

268.01 

280.97 

265.42 

229.70 


177.72 
166.44 
150.93 
171.24 

213.72 
250.67 

261.85 

242.86 
204.51 


150 

201.14 

0.08 

195 

236.92 

240 

246.77 

0.08 

285 

226.71 

mu 

330 

188.72 

0.10 

0 

153.42 

0.10 

15 

143.79 

0.10 

60 

132.00 

0.10 

105 

151.62 

0.10 

150 

191.04 

0.10 

195 

225.33 

0.10 

240 

234.15 

0.10 

285 

213.99 

0.10 

330 

177.361 
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TABLE  D.4:  INSERTION  WINDOW  HYPERSURPACE  (<«8.5*  ,  9.5‘) 


<(#)  •  Ur)  Hp  (Kb)  i(#)  e  fl(#)  Hp  (Km) 


275.29 

265.49 

237.27 

232.56 

261.47 

294.70 

312.63 

310.67 

292.83 


283.08 

268.87 

235.08 


182.33 

171.61 

156.38 

174.43 

214.66 

251.51 

263.36 

245.03 

207.80 


166.91 

157.07 

144.56 

163.05 

202.24 

237.75 

248.15 

228.37 

191.20 


156.75 

147.75 
136.77 
154.86 
192.39 
226.32 
235.50 
215.52 
179.61 
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TABLE  D.5:  INSERTION  WINDOW  HYPERSURFACE  (1-10.5*  .  11.5*) 


<(*)  •  HO  Hp  (Km)  I  |(* )  e  n(*)  Hp  (Km) 


10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 


10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 

10.50 


292.74 

283.92 

255.10 

242.27 

265.30 

298.32 

318.41 

320.37 

307.44 


218.76 

207.91 

187.56 

198.31 

234.60 

271.27 
286.34 
273.67 

242.28 


188.44 

178.47 

164.09 

179.70 

217.20 

253.39 

265.80 

248.22 

212.41 


168.11 

204.72 

239.59 

250.26 

230.94 

194.88 


161.43 

153.06 

143.25 

159.81 

194.91 

228.21 

237.66 

217.77 

182.85 
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TABLE  D.6:  INSERTION  WINDOW  HYPERSURFACE  (i-  12.5*  .  13.5*) 


flC)  «p  (Km) 


»(*) 

e 

cio 

Hp  (Km)  | 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 


12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 


12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 


12.50 

12.50 

12.50 

12.50 


12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 

12.50 


13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 


13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 


13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50! 


13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 


13.50 

13.50 

.3.50 

13.50 

13.50 

13.50 

13.50 

13.50 

13.50 


313.61 

305.67 

276.96 

256.97 
272.16 
303.91 
326.25 
332.82 
325.27 


230.28 

220.39 

200.42 

206.95 

239.11 

274.73 

290.76 

280.01 

251.59 


196.43 

187.31 
173.87 
186.93 
221.15 

256.31 
269.18 
252.45 
218.52 


178.60 

170.41 

160.10 

174.82 

208.59 

242.44 
253.30 

234.44 
199.76 


167.37 
159.81 
151.35 

166.38 
198.87 
231.09 
240.54 
221.01 
187.26 
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«(°) 

e 

no 

Hp  (Km) 

TABLE  D.7:  INSERTION  WINDOW  HYPERSURFACE  .  15.5*) 


{1C)  Hp  (Km) 


337.72 
330.76 

302.73 
277.35 
282.94 
312.04 
336.93 
348.79 

346.73 


244.78 

235.94 

216.55 

218.76 

245.74 

279.82 

296.71 

288.36 

263.40 


206.49 

198.40 
186.18 
196.52 
226.79 
260.54 
273.70 
258.09 

226.41 


186.88 

179.70 

170.87 

183.56 

213.92 

246.40 

257.35 

239.04 


174.93 
168.27 
161.34 
174.66 
204.09 
235.05 
244.50 
225.24 

192.93 


14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 


14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 

14.50 


0.02 

0 

0.02 

15 

0.02 

60 

0.02 

105 

0.02 

150 

0.02 

195 

0.02 

240 

0.02 

285 

0.02 

330 

0.04 

0 

0.04 

15 

0.04 

60 

0.04 

105 

0.04 

150 

0.04 

195 

0.04 

240 

0.04 

285 

0.04 

330 

0.06 

0 

0.06 

15 

0.06 

60 

0.06 

105 

0.06 

150 

0.06 

195 

0.06 

240 

0.06 

285 

330 

0.08 

0 

0.08 

15 

hhbed 

60 

0.08 

105 

0.08 

150 

0.08 

195 

0.08 

240 

285 

0.08 

330 

0.10 

0 

0.10 

15 

0.10 

60 

0.10 

105 

0.10 

150 

0.10 

195 

0.10 

240 

0.10 

285 

0.10 

330 
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Appendix  E 


Supplemental  Figures  for  Section  HI  and  IV 

The  following  pages  contain  3D  graphs  and  contour  mappings  of  the 
tabulated  data  from  Appendix  D. 
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Fig  E.2:  Periapsis  Contour  (i  =  0.5  deg) 


5  deg) 


Fig  E.4:  Periapsis  Contour  (i  =  1.5  deg) 


Fig  E.5:  Per  lope  is  Altitude  (  l  =  2.5  deg) 


Fig  E.6:  Periapsis  Contour  (i  =  2.5  deg) 


Fig  E.7:  Per  lope  le  Al 


100 


Fig  E.8:  Periapsis  Contour  (i  =  3.5  deg) 


Fig  E.9:  Per  lops  is  Altitude  (  l  =  4.5  deg) 


Fig  E.10:  Periopsis  Altitude  (i  =  4.5  deg) 


Is  Altitude  (  l  =  5.5  deg) 


Fig  E.12:  Periapsis  Altitude  (i  =  5.5  deg) 


Fig  E.14:  Periopsis  Altitude  (i  =  6.5  deg) 


F  lg  E.15:  Per  lope  le  Altitude  (  l  =  7.5  deg) 


Fig  E.16:  Periopsis  Altitude  (i  =  7.5  deg) 


p»ni  i  iv 


Fig  E.18:  Periopsis  Altitude  (i  =  8.5  deg) 


.5  deg) 


Fig  E.20:  Periapsis  Altitude  (i  =  9.5  deg) 


F  lg  E. 21  :  Per  lo 


Fig  E.22:  Periapsis  Altitude  (i  =  10.5  deg) 


Fig  E.23:  Per  lops  Is  Altitude  (  l  =  11.5  deg) 


24:  Periapsis  Altitude  (i  =  11.5  deg) 


Fig  E.26:  Periapsis  Altitude  (i  =  12.5  deg) 


Fig  E.27:  Per  lops  le  Altitude  (  l  =  13.5  deg) 


Fig  E.28:  Periapsis  Altitude  (i  =  13.5  deg) 


Fig  E.30:  Periapsis  Attitude  (i  =  14.5  deg) 


Fig  E.32:  Periopsis  Altitude  (i  =  15.5  deg) 


Fig  E.34:  Periapsis  Contour  (e  =  0.02) 


Altitude  Ce  =  0.04) 


Fig  E.36:  Periapsis  Contour  (e  =  0.04) 
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Fig  E.38:  Periapsis  Contour  (e  =  0.06) 


(wjf  )  mpni)t/y  •/  m&o) 


Fig  E.40:  Periopsis  Contour  (e  =  0.08) 


Flo  E.41:  Per  lope  le  Altitude  (e  =  0.10) 


Fig  E.42:  Periapsis  Contour  (e  =  0.10) 


Fig  E.43:  Par  laps  l®  Altitude  (  =0  deg) 


Fig  E.44:  Periapsis  Contour  (  =0  deg) 


F ig  E. 45:  P®r  lo 


Fig  E.46:  Periapsis  Contour  (  =15  deg) 


Fig  E.47:  Per  lops  Is  Altitude  (  =  60  deg) 


Fig  E.48:  Periapsis  Contour  (  =60  deg) 


Fig  E.  49:  Per  lops  Is  Altitude  (  =  105  deg) 


Fig  E.50:  Periapsis  Contour  (  =105  deg) 


F ig  E. 51  :  Per  lo 


Fig  E.52:  Periapsis  Contour  (  =150  deg) 


Fig  E. 53:  Por  lop 


Fig  E.54:  Periapsis  Contour  (  =195  deg) 


Fig  E.55:  Per  lope  Is  Altitude  C  =  240  deg) 


Fig  E.56:  Periapsis  Contour  (  =  240  deg) 


Fig  E.57:  Por  lo 


Fig  E.58:  Periapsis  Contour  (  =  285  deg) 


Fig  E.60:  Periapsis  Contour  (  =  330  deg) 
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